1. Field of the Invention
This invention relates to digital filters.
2. Description of the Prior Art
In digital signal processing apparatus such as a digital audio mixing console, recursive digital filters are used for signal processing operations such as those required in equalisation sections.
The general form of a transfer response of a first order recursive digital filter is:                               H          ⁡                      (            z            )                          =                                            a              0                        +                                          a                1                            ⁢                              z                                  -                  1                                                                          1            +                                          b                1                            ⁢                              z                                  -                  1                                                                                        (        1        )            
This equation is implemented digitally as a series of multiplies and delays, corresponding to each term in the equation.
An example implementation of the filter of equation (1) is shown in FIG. 1 of the accompanying drawings. In FIG. 1, the input signal samples are supplied to a multiplier 10 (multiplication by a coefficient a0) and also to a one-sample delay (zxe2x88x921) 20. Samples delayed by the delay 20 are multiplied 30 by a coefficient a1. An adder 40 adds the output of the multipliers 10, 30 and the output of a further multiplier 60 which multiplies output samples y(n), delayed 50 by one sample, by a coefficient b1. The output of the adder 40 forms the output samples y(n).
A problem with the use of such filters is that the so-called limit cycle behaviour can reduce the performance of the filters in high fidelity audio applications.
One example of this limit cycle behaviour is the response of a filter when the input of the filter is reduced to zero. Rounding errors caused by the use of finite precision arithmetic mean that the output of the filter does not necessarily reduce to zero when the input reduces to zero. Instead, the output can stay at a non-zero value or oscillate about zero. In an audio system this is heard as a deterioration of the audio performance.
An example of limit cycle behaviour can be derived from equation 1. The difference equation by which output samples y(n) may be calculated from input samples x(n) is as follows:
y(n)=a0x(n)+a1x(nxe2x88x921)xe2x88x92b1y(nxe2x88x921)xe2x80x83xe2x80x83(2)
To demonstrate the zero input limit cycle behaviour, x(n) and x(nxe2x88x921) are set to zero and an initial condition is substituted for y(nxe2x88x921):
y(n)=0xe2x88x92b1y(nxe2x88x921)xe2x80x83xe2x80x83(3)
In a specific example, for a filter having a sampling rate of 25Hz and a decay time of approximately 2 seconds, a0=a1=0.045, and b1=xe2x88x920.955. Taking the initial condition of y(xe2x88x921)=15, the following outputs for integer arithmetic can be derived:
This example demonstrates a limit cycle with zero frequency (i.e. a steady error value is output). However, if b1 were positive then the rounded y(n) would alternate in sign, oscillating at one half of the sampling rate. As the absolute value of b1 approaches 1, the magnitude of the zero input limit cycle output gets larger and larger.
Attempts to overcome or alleviate this problem have included varying the circuit or program design used to implement the filterxe2x80x94for example, the so-called xe2x80x9cLattice formxe2x80x9d, the xe2x80x9cTranspose formxe2x80x9d, the xe2x80x9cGold-Rader Structurexe2x80x9d and others have been tried, but even where these filters have given a slightly improved limit cycle performance, they have had a worse performance in another area important in audio applications.
Double precision operation has also been tried, as a straightforward way of reducing the effect of the rounding errors which lead to limit cycle problems, but their use requires too much processing capacity to be a useful solution. In some instances, where double precision operations (e.g. multiplies or additions) are built from combinations of single precision operations, it has been found that a double precision filter requires seven times the number of operations of an equivalent single precision filter.
This invention provides a recursive nth-order digital filter, where n is greater than 1, the filter having a pole-zero pair substantially co-located in the z-plane.
The elegantly simple solution provided by the invention is to change the way in which a filter response is calculated (i.e. implemented as a digital filter) without necessarily changing the actual response of the filter.
Referring to the example of equation (1) above, the response described in equation (1) is multiplied by a further expression having an identical numerator and denominator (i.e. a further expression equal to one):                               H          ⁡                      (            z            )                          =                                                            a                0                            +                                                a                  1                                ⁢                                  z                                      -                    1                                                                                      1              +                                                b                  1                                ⁢                                  z                                      -                    1                                                                                ·                                    1              +                              γ                ⁢                                  xe2x80x83                                ⁢                                  z                                      -                    1                                                                                      1              +                              γ                ⁢                                  xe2x80x83                                ⁢                                  z                                      -                    1                                                                                                          (        4        )            
This gives a revised expression which is equal to the first equation given above.                               H          ⁡                      (            z            )                          =                                            a              0                        +                                          z                                  -                  1                                            ⁡                              (                                                                            a                      0                                        ⁢                    γ                                    +                                      a                    1                                                  )                                      +                          γ              ⁢                              xe2x80x83                            ⁢                              a                1                            ⁢                              z                                  -                  2                                                                          1            +                                          z                                  -                  1                                            ⁡                              (                                                      b                    1                                    +                  γ                                )                                      +                                          b                1                            ⁢              γ              ⁢                              xe2x80x83                            ⁢                              z                                  -                  2                                                                                        (        5        )            
Because this expression is equal to the first expression given above, it provides the same filter response. However, more terms are used in calculating the response, which has the effect that when it is implemented digitally, limit cycle problems caused by rounding errors are reduced.
Another view of the revised filter implementation is that the equivalent of an additional pole and an additional zero have been added at identical positions in the z plane. The filter of equation (5) is now a second order filter but has the same frequency and phase response (except in limit cycle behaviour) as the first order filter of equation (1). The additional processing overhead required to implement the revised filter of equation (5) is much less than the seven-fold increase needed to use double precision arithmetic. In embodiments of the invention (e.g. an audio filter embodiment to be described below) an improvement of around 10 dB in the signal to noise ratio can be obtained at audio frequencies below about 100 Hz.
The invention also provides a method of generating filter coefficients to provide a desired nth order filter response (e.g. where n is an integer greater than or equal to 1), the method comprising adding one or more substantially co-located pole and zero pairs to form an (n+m)th order filter having the desired nth order filter response (e.g. where m is an integer equal to 1 or more).